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Abstract 

During direct numerical simulation of the isotropic turbulence of surface gravity 
waves in the framework of Hamiltonian equations formation of the long wave back- 
ground or condensate was observed. Exponents of the direct cascade spectra at the 
different levels of an artificial condensate suppression show a tendency to become 
closer to the prediction of the wave turbulence theory at lower levels of condensate. 
A simple qualitative explanation of the mechanism of this phenomenon is proposed. 
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1 Introduction 



Description of the water waves appeals scientists attention during several cen- 
turies. At the same time, first attempts to explain observed spatial and tem- 
poral spectra are relatively recent. One of the first and at the same time most 
famous works is the paper by Phillips pQ which was, probably, the first attempt 
to give an explanation for power-like spectra of surface gravity waves observed 
in numerous experiments. In just one decade statistical theory of water waves, 
based on the the kinetic equation for waves derived by Hasselmann [2J and 
solutions of this equation obtained from Zakharov's theory of wave (or weak) 
turbulence [3114"] . These solutions [5.4J are stationary Kolmogorov solutions 
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of the kinetic equation corresponding to flux of energy from large to small 
scales (direct cascade) and flux of wave action (waves "number") from small 
to large scales (inverse cascade). Now the kinetic equation is a base tool for 
wave forecasting. 

The Hasselmann kinetic equation for waves was derived under some assump- 
tions, which include Gaussian statistics for the wave field and prevalence of 
resonant interactions [3]. These assumptions are subject to confirmation. Mod- 
ern numerical methods allow us to perform wave field modeling in the frame- 
work of kinetic equations faster than the real processes in nature. At the same 
time, it is still not possible to create a wave forecasting model based on di- 
rect numerical simulation of the primordial dynamic equations. Fortunately, 
in practical applications we do not need to know velocity and elevation at ev- 
ery point of the surface. Statistics, especially mean wave hight and speed, are 
what really matters for estimation of operational conditions of oil platforms 
and cargo ships. Such statistics are exactly the subject of the theory of weak 
turbulence. As a consequence, the problem of confirmation and correction of 
the theory of wave turbulence is of great practical importance. 

Open sea observations confirm temporal and space spectra predicted by the 
theory of wave turbulence [6.7.8J. A comprehensive review of experiments and 
comparison with the theory of the weak turbulence can be found in [9lll0fllj . 
At the same time even the most recent, state of the art wave tanks experiments 
give strange results [T2~|IT5] . demonstrating dependence of the spectrum slope 
on the average steepness or forcing level. 

One of the ways from this deadend is direct numerical simulation in the frame- 
work of dynamical equations. Such numerical experiments can be separated 
in two major groups. In swell simulation experiments initial condition is some 
spectral distribution, there is no pumping. Temporal evolution of such initial 
condition is simulated and than can be compared with corresponding simula- 
tion in the framework of kinetic equation |14lll5fl6fl7lll8fl9] . More interesting 
and more complex is the simulation of the wave turbulence with pumping. 
Classical sandbox for such simulations is an isotropic turbulence, which al- 
lows us to demonstrate of the features of the wave turbulence and, at the 
same time, simplifies situation due to additional smoothing of the spectra af- 
ter angle averaging. Although this situation looks pretty artificial, the isotropic 
turbulence can be observed in the sea with large amount of broken ice, for ex- 
ample, where waves create isotropic wavefield after multiple reflections. But in 
such a case floating ice have to be taken into account, which makes situation 
more complex. Simulation of isotropic wave turbulence already confirmed some 
prediction of the theory of weak turbulence. Direct cascade spectra were ob- 
tained in [20|21f22| . Inverse cascade spectra, although in simplified equations, 
were obtained in [23]. For the first time simultaneous simulation of direct and 
inverse cascades for gravity waves in the framework of primordial dynamical 
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equations was performed in [2 



Current paper is a continuation of the work [21]. The present article will 
concentrate on the condensate generation and its influence on the direct cas- 
cade spectra of wave turbulence of surface gravity waves. A more detailed 
and deep insight in the mechanism of the condensate influence is presented. 
Starting from the results of isotropic turbulence simulation author will explore 
in details how presence of the condensate changes the slope of direct cascade 
spetrum. It will be demostrated that the presence of the inverse cascade is also 
important. At the same time in the pure situation with suppressed conden- 
sate and inverse cascade observation shows spectrum slope close to previous 
numerous results of simulations [2"U|21ll22j . A simple qualitative explanation 
of the mechanism of spectrum distortion is provided. Direct analysis of the 
simulation results supports the theory. 



2 Theoretical background. 

Author follows previous works in theoretical description of the system under 
consideration l20.2TlfT9l[2H. 



2.1 Dynamical equations. 

We study the potential flow of an ideal inviscid incompressible fluid with the 
velocity potential (p = <p(x, y, z; t) 



in the infinitely deep domain occupied by the fluid. Equations for the boundary 
conditions at the surface are the following 



Here i](x,y;t) is the surface elevation with respect to still water, g is the 
gravity acceleration. Equations (pQ) are Hamiltonian [25] with the canonical 
variables r)(x, y; t) and ip(x, y\ t) = <f>(x, y, rj(x, y\ t); t) 



A<p = 0, 




(1) 



dr] 5H dip 
~di~lnp' ~dt 



5H 



8r] ' 



(2) 
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where H is the Hamiltonian of the system 

+00 / T) 

H = - J dxdylg V 2 + J \V<P\ 2 dz 

—00 \ —00 j 

Unfortunately H cannot be written in the close form as a functional of rj and 
ip. However one can limit Hamiltonian by first three terms of expansion in 
powers of 77 and tp 



H = H 
Hi = - 



grf + ipkip ) dxdy, 



rj 



\Vip\ 2 - {kip) 2 dxdy 



(3) 



H2 = - j f)(kip) k(r](kip)) + TjAip dxdy. 



Here k is a linear integral operator (k = ^—Aj , such that in fc-space it cor- 
responds to multiplication of Fourier harmonics (ip^ = ^ J ip?e lkr dxdy) by 
k — \/k 2 + k%. For gravity waves this reduced Hamiltonian describes four- 



acquire the form 



wave interaction. 

Then dynamical equations 



fj = kip — (V(r]'Vip)) — k[r]kip] + 
+k{ v k[ v kip}) + ^A[rfkip\ + h[ v 2 AiP} - F-^rik], 

i> = -m-\ [(W) 2 - (h) 2 } - (4) 

-[kiP]k[ V kiP] - [ V ktP]AtP - F-^iPk] + F -1 [Pg]. 

Here dot means time-derivative, A — Laplace operator, A; is a linear integral 
operator (jfe = ^A), F~ l is an inverse Fourier transform, 7^ is a dissipation 
rate (according to recent work |26j it has to be included in both equations), 
which corresponds to viscosity on small scales and, if needed, "artificial" damp- 
ing on large scales. P^ is the driving term which simulates pumping on large 
scales (for example, due to wind). In the fc-space supports of 7^ and P^ are 
separated by the inertial interval, where the Kolmogorov-type solution can be 
recognized. These equations were derived as a result of Hamiltonian expansion 
in terms of krj. From the physical point of view the fc-operator is close to the 
derivative operator, so we expand in powers of the slope of the surface. In 
most of the experimental observations average slope of the open sea surface /1 
is of the order of 0.1, so such expansion is very reasonable. 
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These are dynamical equations which will be subject for modeling in a signif- 
icant part of this article. In order to understand the motivation of the paper 
we need to introduce statistical description of the wave field. 



2.2 Kinetic equation. 



For the analisis of the wave field in is natural to use Fourier representation: 



Wave field functions i/)(r,t) and r)(r,t) are real valued functions, so their 
Fourier transforms have Hermitian symmetry 



As a result we have two complex functions and each of them use effectively 
only half of the wavenumbers plane. For further derivation is is convenient to 
introduce so called normal variables: 



In this variables for a wave field of one plain wave with wavevector ko only 
harmonics a^ Q will be nonzero. We can state that we introduced elementary ex- 
citations of the system. Hamiltonian equations in this variables has canonical 
form: 



^ = 2^ / Vve ifc> d 2 r, Vi=—J Vr ^d 2 r. 




(5) 



k k 



(6) 
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Hamiltonian in this terms allows us to separate different wave interaction 
processes: 



H = LU k \a^\ 2 dk, 

lift - - - - 

Hl = 62^J E ikS a ki a k2 a k Q + a i a h a W 5 ( kl + k2 + k o)dhdk 2 dk + 

\ \ f i* -* —> 

2 2^ J C UM^°% + a k a h a k ^ kl +k2 ~ ko)dhdk 2 dk , 
4(2^ / W ^^ a ^ a ^ + a k a h a k a V x 



(7) 



1 1 

+ ; 

H 2 

xS(ki + k 2 + k 3 + k4)dkidk 2 dk 3 dk 4: + 
1 1 r 

+—-, — — / F? t r r (a>r clt &t a>z + «r o*r a*? a?)x 

4 (2n) 2 J kik 2 k 3 k4\ ki k 2 k 3 ki fci k 2 k z ki' 

xS(ki — k 2 — k 3 — fc^d^idA^d/i^d/^-l- 

1 1 r ^^^^^^^^ 
+ 4?>)2 J D k^kS a ki a k2 a l 3 a *k4 6 ^ kl + k 2- k 3- h)dhdk 2 dk 3 dk 4 . 



We omit formulae for matrix elements of these processes (they can be found 
in [27]) for the sake of conciseness. 

For gravity waves on the surface of the deep fluid dispersion relation is u = 
y/gk and resonant conditions 



u ko = u kl + u k , 2 , k = ki + k 2 , 

are never fulfilled, which means that three-waves processes of decaying and 
merging of the waves are prohibited. Such dispersion relations are of the "non- 
decay type" . 

In this case it is well known fact that cubic terms of the Hamiltonian flTj) can be 
eliminated by proper canonical transformation a(k, t) — > b(k, t). Hamiltonian 
in new variables looks as follows 



H = J uj k \b^\ 2 dk 



fci 

#i = 0, 

1 1 r ^^^^^^^^ 

#2 = J Tk^ksk^h^kshA^ + k 2~h- h)dhdk 2 dk 3 dh. 

The formula of the matrix element T^^^ is too bulky to give it here and 
can be found in |28|. 
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If we suppose, that in our field phases and amplitudes are not correlated one 
can introduce pair correlator 

(a^,) = n k 5(k-k'). (8) 



The ri£ is measurable quantity, connected directly with observable correlation 
functions. For instance, from definition one can get 

4HM 2 > = ^K + n_ fc ). (9) 
z 9 



For gravity waves we can also use the following correlation function 

(b n bt) = N k 6(k-k'). (10) 

Because in the case of variables all non-resonant processes are excluded this 
is the correlation function for which kinetic equation is derived. 

The relation connecting n k and Nk is very simple (in the case of deep wa- 
ter) [29] 

n k ~N k 

11 

n k 

where \x = (ka) 2 , here a is a characteristic elevation of the free surface. In the 
case of the weak turbulence \i <C 1. 

The correlation function N k obeys the kinetic equation [2"t5fl] 
dN k 



dt 



S nl (N, N, N) + S in (k) - S diss (k), (12) 



Here Si n — corresponds to input of energy, Sdi SS — dissipation of energy in 
the system, and 



7^- .. T 



2 

X 



S nl (N, N J N)=4ir 
x(N kl N k2 N k3 + N k N k2 N k3 - N k N kl N k2 - (13) 
-N k N kl N k3 )5(k + ki — k 2 — k 3 )5(u k + u kl - uj k2 - wjjdfcidfodfcs, 



describes nonlinear interaction of waves. Expression for £ % ^ can be found 
in EH. 



7 



2.3 Kolmogorov-Zakharov solutions of kinetic equation. 

Let us consider stationary solutions of the kinetic equation assuming that 

• The medium is invariant with respect to rotations; 

• Dispersion relation is a power-like function uo = ak a ; 

• T kMMM is a homogeneous function: T ^ ^ ^ = e^T^^. 

Under these assumptions Zakharov [5, 3,30,4] obtained Kolmogorov-like solu- 
tions corresponding to fluxes of two intergals of motion (energy and wave 
action or number of waves): 

4 2) = C 2 Q^k- 2J ^~ d 

Here d is a spatial dimension (d = 2 in our case). In the case of gravity waves 
on a deep water u = \fgk (a = 1/2) and, apparently, (3 = 3. As a result one 
can get: 

4 2) = C 2 Q 1 /3 A .-23/6 

The first spectrum njp corresponds to the direct cascade, describing flux of 
energy from large scales (pumping) to small scales (dissipation). The second 
spectrum n k describes to inverse cascade, corresponding to flux of number of 
waves (or wave action) from small scales (pumping) to larger scales. In the case 
of swell or decaying turbulence inverse cascade reveals itself as a downshift of 
spectral maximum to smaller wavenumbers, resulting in growing long waves. 

2.4 Verification of the kinetic equation. 

As it was shown during derivation of kinetic equation, we neglect some pro- 
cesses, such as multiple harmonics generation, because they are nonresonant. 
Also, statistical description by definition does not take into account rare or 
catastrophic events, for example solitons, freak wave formation, and wave 
breaking or whitecapping which could be very important. It means that we 
need to verify these assumptions using real life or numerical experiments. 

Broad program of verification of kinetic equation for gravity waves was fulfilled 
during the last decade. The key method of verification is numerical simula- 
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(14) 



(15) 



tion in the framework of primordial Hamiltonian equations for shape of the 
surface T](f,t) and velocity potential on the surface if)(r,t). Such numerical 
experiments provide complete information about the system, including phase 
and amplitude of calculated function at every points of the surface at every 
moment of time. Of course such luxury comes at a great cost of enormous 
computational complexity. Due to rapid progress of computers available for 
computations and numerical algorithms during last couple of decades such 
simulations became possible. The first work which presented successful simu- 
lation of weak turbulence of gravity waves was published in the beginning of 
2000s [T5] . It was simulation of the decaying turbulence. Decaying turbulence 
experiments start from initial spectral distribution close to realistic swell which 
is propagating without pumping, decaying (loosing energy) due to nonlinear 
transfer from spectral peak to the dissipative scale. In the paper [15] forma- 
tion of weak turbulent Kolmogorov-Zakharov [4] spectral tail was confirmed. 
In the series of works [2"0l27f21] the first confirmation of Kolmogorov-Zakharov 
spectra in systems with pumping was provided as well as confirmation of in- 
dependance of spectrum exponent from the pumping force in some reasonable 
interval of input fluxes. Shortly after that some possible obstacles for verifi- 
cation of Hasselmann equation through wave tank experiments were pointed 
out [TTJ . In papers [T9"1[TS] idea of modification of the dissipation term Sdi SS 
in order to take additional energy transfer to the dissipation scales due to 
whitecapping and wavebreaking was developed through a series of massive 
numerical simulations of decaying turbulence. In a recent paper [23] some 
consequences of whitecapping influence on a Kolmogorov-Zakharov spectra in 
wavetanks experiments were analyzed. 



3 Numerical simulation 

3. 1 Simulation parameters. 

Here we continue investigation of the problem which was started in [24]. We 
simulated the primordial dynamical equations (T4|) in a periodic spatial domain 
27r x 27T. In order to save computation time multigrid approach has been 
used. Because formation of the inverse cascade is the most time consuming 
part of calculations (characteristic nonlinear time is about Uk/fi — lOu^ and 
&k ~ vk), we completed this calculation on the relatively small grid 256 x 256. 
After that we gradually increased resolution up to 1024 x 1024 nodes. To check 
whether we loose something in the high frequency region long time simulation 
on grid 256 x 256 was performed. It demostrated absence of change in the 
pumping and high frequencies region. 

The used numerical code was verified in |31II20|21|I17|I18|19] . Gravity acceler- 
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ation was g = 1. The pseudo- viscous damping coefficient had the following 
form 

Ho°^l (16) 



where k d = 256 and 70,1024 = 2.7 x 10 4 for the 1024 x 1024 grid and k d = 64 
and 70,256 = 2.4 x 10 2 for the smaller 256 x 256 grid. Pumping was an isotropic 
driving force narrow in wavenumber space with random phase: 

{> 77, (fc-fcpi)(fc p 2-fc) 
fe 2 -V) 2 ' (17) 
— if k < k p i or k > /c p2 ; 



here k p \ = 28, k P 2 = 32 and Fq = 1.5 x 10 5 ; R%(t) was a uniformly distributed 
random number in the interval (0, 27r] for each k and t. The initial condition 
was low amplitude noise in all harmonics. Time steps were Atio24 = 6.7 x 10~ 4 
and At256 = 5.0 x 10 ~ 3 . We used Fourier series in the following form: 



2tt 2tt 

1 



% = = —2j J Vr^d\ 





(2tt 



N x /2-l Ny/2-1 
lr m 1 \ \ * ifcr 



f- 1 IvA= E E vn 



c 



-TV^/2 



As it was reported in [21j as a result of simulation we observed formation of 
both direct and inverse cascades (Fig. [TJ solid line), although exponents of 
power-like spectra were different from weak turbulent solutions (fl"5l) . It is im- 
portant to note that development of the inverse cascade spectrum was arrested 
by the discreteness of the wavenumber grid in agreement with [3"T1I17II32|I3"3"] . 
After that, a large scale condensate started to form. The reason of this phe- 
nomenon is the following. It can be shown that resonant conditions for four 
waves interaction of gravity waves are never fulfilled exactly on the discrete 
grid of wavevectors (integer number for 27r x 27r box): 

^fci + wjfc 2 ~ ^3 ~ ^4 7^ 0- 

Nonlinear frequency shift T weakens resonant condition, now resonant curve 
has finite width: 

Uki + UJ k 2 ~ ^hj, — < r. 
Nonlinear frequency shift depends on Tjg^g j^, which is homogeneous func- 
tion. It grows as k 3 when k is increased. This is the reason why it is much easier 
to get resonant interaction in the direct cascade region (high wavenumbers) . It 
also decreases as k 3 when k — > 0. It means that effectively grid become more 
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and more coarse at large scales. It means that at some scale inverse cascade is 
stopped because four-waves nonlinear interaction is "turned off" . At the same 
time flux of action still brings new waves to this scale. We have "condensa- 
tion" [31] of waves. Similar processes were observed in 2D-hydrodynamics [35] 
and plasma [36] 

As one can see, the value of wave action |afc| 2 at the condensate region is more 
than an order of magnitude higher than for the closest harmonic of the in- 
verse cascade. The dynamics of large scales became extremely slow after this 
point. We managed to achieve downshift of the condensate peak by one step of 
the wavenumber grid during long time simulation on a small grid (256 x 256) 
(Fig. [TJ line with long dashes). As one can see, we observed elongation of the 
inverse cascade interval without significant change of the slope. Unfortunately, 
the inertial interval for inverse cascade was too short to exclude the possible 
influence of pumping and condensate. For direct cascade spectra we also used 
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Fig. 1. Spectra (|afc| 2 ). With condensate on the 1024 x 1024 grid (solid); on the 
256 x 256 grid with more developed condensate (dashed); without condensate on 
the 1024 x 1024 grid (dotted); without inverse cascade on the 1024 x 1024 grid 
(dashed-dotted). 

a log- log scale. Results are present in Fig. [3j (left). Formally, in this case we 
have quite a long inertial interval 32 < k < 256, but in reality damping has an 
influence on the spectrum approximately up to k ~ 180. Still in this case we 
have more than half of a decade. The theory of weak turbulence gives us depen- 
dence ~ k~ 4 (jl5]) . known as Kolmogorov-Zakharov direct cascade spectrum 
for gravity waves. Nevertheless, one can see that we observe k~ 9 ^ 2 , known 
as Phillips [TTI3T] spectrum. Possible explanation was provided in [21]. The 



11 



main idea was the following: presence of the condensate stimulated additional 
whitecapping and wavebraking processes by local increase of steepness. These 
processes provide additional transfer of energy to the dissipation region. When 
this dissipation is strong enough to balance flux of energy, such mechanism 
should give Phillips spectrum as a result [37]. It was shown, that suppression 
of the condensate (see Figure [2]) changes spectrum exponent to some value 





Fig. 2. Surface of the fluid rj(r) with (left) and without (right) condensate. 



closer to the theory of weak turbulence. However, at that time we failed to 
avoid long waves' influence on the direct cascade spectrum completely. The 
longest waves of inverse cascade were providing long wave background which 
influence distorted the spectrum. 



In order to eliminate this influence let us suppress inverse cascade completely 
by artificial damping in low wavenumbers. The result of simulation is shown 
in Figure [3] (right panel) . 




Fig. 3. Compensated direct cascade spectra C(\ak\ 2 )k s with condensate and inverse 
cascade (left) and without them (right). 



For comparison we have put in the same figure initial spectrum for direct 
cascade (left). As one can see, the exponent of the spectrum changed and now 
almost coincides with the prediction of the weak turbulent theory. The tiny 
difference may be a result of the influence of the left peak, which is the result 
of insufficient length of damping region. At the same time, longer artificial 
damping can start to influence the pumping area and change the input flux 
significantly. 
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backgroung periodic current 



Fig. 4. Initial wave (solid line) become steeper (dashed line) in the result of inter- 
action with the backround current. 

A qualitative explanation of the condensate's influence on the short waves 
could be the following: let us consider a propagating wave with some given 
slope at its front; a much longer wave can be treated as a presence of a strong 
background flow. If the direction of the flow is opposite to direction of the 
wave's propagation, the slope of the wave's front will increase (see Figure HJ). 

This is what we see in our simulations. The average steepness \i = \J (|Vr?| 2 ) 
has changed: with condensate \i c ~ 0.14 (Figure EJ), without condensate fi nc ~ 
0.12 (Figure [6]), after suppression of inverse cascade \i nc ~ 0.1 (Figure [7]). We 
should mention that very high values of gradient of the surface (higher that 
value for a limiting Stokes wave) are possible in our model only due to very 
strong damping at high wavenumbers and relatively narrow inertial interval, 
which prevents waves from real breaking by strong dissipation of energy for 
every breaking onset. We regularize our model by this dissipation. At the same 
time here we don't need to know details of the energy dissipation (how in 
details wave breaks and produce very short waves), all we need is to dissipate 
all the energy involved in the breaking. 



4 Conclusion 



The observed phenomena of condensate formation and its influence on the 
direct cascade have to be quite common for experiments in experimental wave 
tanks. Finite size of the surface area gives us the same discrete grid of wavevec- 
tors as in the case of periodic boundary conditons. It means that at some 
large scale four waves nonlinear interaction will be arrested, which provides 
conditions for condensate formation. Such condensate will distort the direct 
cascade spectrum, changing its exponent through stimulation of whitecapping, 
resulting in additional dissipation. It is worth to note that this is threshold 
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Fig. 5. Gradient of some part of the surface. In the presence of condensate. 




Fig. 6. Gradient of some part of the surface. In the absence of condensate. 

phenomenon with respect to average steepness as it was shown in [38]. If 
such dissipation is strong enough to balance the flux of energy, following [37] 
we shall get Phillips spectrum, as it was demonstrated in the present paper. 
Gradually suppressing condensate and inverse cascade one can achieve spec- 
trum corresponding to the stationary solution of the kinetic equation. Simple 
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Y 



Fig. 7. Gradient of some part of the surface. In the absence of inverse cascade. 



mechanism of condensate interaction with short waves of direct cascade was 
proposed. This interaction along with other distortion possibilities (see, for 
example, [39]) can be the reason of observed spectra deviation from the weak 
turbulence predictions. These results can be used for explanation of experi- 
mental data obtained in the wavetanks. 
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